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A model is presented for the contrarotating vortex pair that is 
formed by a round, turbulent, subsonic jet directed normally into a 
uniform, subsonic cross flow. The model consists of a set of algebraic 
equations that describe the properties of the vortex pair as a function 
of their location in the jet plume. The parameters of the model are 
physical characteristics of the vortices such as the vortex strength, 
spacing and core size. These parameters are determined by velocity 
measurements at selective points in the jet plume. 
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INTRODUCTION 


interesting aerodynamic problem occurs when a vertical takeoff 
and landing aircraft transitions from hovering to forward flight. The 
cross flow caused by the aircraft's forward flight interacts with the 
aircraft’s lifting jets, and brings about a loss of performance in 
addition to stability problems. Figure 1 is a sketch of a VTOL air- 
craft with a single lifting jet transitioning from hovering to forward 
flight. The interference between the lifting jet and the cross flow 
induces a low pressure region in the jet wake and a distribution of 
pressure induced forces over the aircraft as indicated schematically 
in the figure. The flow field associated with a VTOL aircraft and its 
multiple lifting jets is very complex. In order to simplify the 
problem, while still retaining the essential characteristics of the 
jet and cross flow interaction process, it has been customary to 
restrict the problem to that of a single, subsonic jet exhausting 
normally through a flat plate into a uniform, subsonic cross flow. 

Other applications of the jet in a cross flow are the cooling of jet 
turbine combustors, where coolant gases are injected into the combustors 
to dissipate heat, and the environmental problem of the discharge of 
cooling water from a power plant into a waterway. 

Early investigations into the jet in a cross flow were concerned 
with a qualitative description of the flow. Several of these investi- 
gations resulted in empirical equations for the jet centerline, which 
was defined to be the locus of maximum velocity in the symmetry plane. 
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Many of the flow characteristics, it was found, depended on the jet to 
cross flow momentum ratio. Several investigators found it convenient 
to use the square root of the momentum ratio and defined it as the 
effective velocity ratio. In certain restrictive cases, it was shown 
that the effective velocity ratio reduced to the jet to cross flow 
velocity ratio, (ref. 1). 

Recently there have been several detailed studies into the veloc- 
ity and pressure fields associated with the jet in a cross flow. In 
each of these investigations it was found that the jet wake was dominated 
by a pair of contrarotating vortices. Kamotani and Greber (ref. 2) and 
Harms (ref. 3) investigated the effect of heated jets on the flow field. 
Thompson (ref. 4), utilizing an unheated jet, conducted a study of the 
pressure distribution on a flat plate and the velocity field induced by 
a jet in a cross flow. Fearn and Weston (ref. 1 and 5) investigated the 
velocity and pressure fields associated with an unheated jet directed 
normally into a cross flow. Fearn and Weston presented simple analytic 
models for the vortex pair associated with a jet in a cross flow. 

This paper will investigate the effect of varying certain parameters 
in one of the models presented by Fearn and Weston. The vortex model 
will also be extended to give a more convenient means of determining 


the vortex properties. 



LITERATURE REVIEW 


An important feature of the jet in a cross flow is the contra- 
rotating vortices that roll up in the wake region. Figure 2 is a 
sketch of the wake region induced by a jet in a cross flow. A plane 
of symmetry (Y=0) is seen to exist in the flow field. The vortex 
centers are defined as the points of maximum vorticity in a cross 
section and are shown schematically in figure 2. The loci of the 
vortex centers are commonly referred to as the vortex trajectories. 

It is convenient to define the projections of the vortex trajectories 
onto the symmetry plane as the vortex curve, shown schematically in 
figure 3. The vortex curve is then described as the locus of the 
midpoints of the line joining the vortex centers. It is seen in 
figure 3 that the vortex curve falls slightly below the jet centerline. 
Figure 3 also illustrates the reference frames commonly in use. One 
reference frame has its origin at the center of the jet exit and is 
aligned with the tunnel coordinate system. Two other coordinate sys- 
tems are aligned with cross sections perpendicular to the jet center- 
line and the vortex curve. The system referred to the jet certerline 
has its axes denoted with a subscript j while the system referred to 
the vortex curve has its axes denoted with a subscript v. 

Kamotani and Greber (ref. 2), who were prompted by the gas turbine 
combustor cooling problem, utilized a 1/4 inch diameter jet in their 
study of heated and unheated jets in a cross flow. The authors found 
that the flow was dominated by a pair of contrarotating vortices that 
form behind the jet and persist for a long distance downstream of the 
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jet orifice. This compares with the findings of Pratte and Baines 
(ref. 6), in which the vortices were located as far as 1000 jet diam- 
eters downstream of the jet orifice. Kamotani and Greber presented 
the data as plots of velocity and temperature distributions in the 
symmetry plane and- in planes perpendicular to the jet centerline. The 
vortices were apparent when the velocity vectors in the wake region 
were projected into cross section planes perpendicular to the jet 
centerline. Results were shown for effective velocity ratios of 3.91 
and 7.72. The authors made no attempt to calculate the strength of the 
vortex pair. 

The aerodynamic problems associated with V/STOL aircraft have 
motivated several recent papers. Harms (ref. 3) studied the temperature 
effects on the velocity field induced by a 5 centimeter diameter jet 
issuing into a cross flow. Extensive velocity measurements were taken 
in planes perpendicular to the tunnel axis for an effective velocity 
ratio of 8. It was found that the position of the vortex centers re- 
mained essentially the same for hot and cold jets of the same effective 
velocity ratio. The only difference being that for the hot jet the 
vortices were more diffuse. Harms stated that the vortex pair absorbed 
the axial momentum of the jet and were dissipated in the far field by 
the action of viscous forces. Harms’ study was similar to Kamotani and 
Greber in that he recognized that the vortices were the dominant feature, 
but no attempt was made to calculate their strength. 

Thompson (ref. 4) investigated the ground board pressure distribu- 
tion and velocity field induced by unheated jets issuing into a cross 
flow. Thompson utilized a one inch diameter jet, in addition to elliptic 
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jets of comparable exit area. Velocity measurements were taken in 
cross section planes perpendicular to the jet centerline. 

Thompson attempted to infer the strength and location of the 
vortices induced by a circular jet for effective velocity ratios of 
2, 4 and 8. The vortex centers were located from the distribution 
of the sidewash component (V ) , of the velocity in cross section 
planes approximately perpendicular to the vortex curve. The dis- 
tributions were found by taking velocity measurements in traverses 
parallel to the and Y^ axes in the cross section plane. Figure 4a 
is a sketch illustrating typical distributions in the cross section 
planes. The component changed sign as the line joining the vortex 
centers was traversed in a Z direction. The vortex separation was 
determined by locating the peaks in the distribution as the in- 
dividual vortices were traversed in a direction. Thompson assumed 
that the vortex properties changed slowly enough with that the 
vortices could be treated as though they were two-dimensional. The 
velocity induced by a single, two-dimensional vortex filament was equal 
to r/2-rrr , where F was the vortex strength and r was the. radial distance 
from the vortex center. Thompson measured the component of velocity 
in a traverse parallel to the Z^ axis in the cross section and passing 
through the vortex center. Figure 4b is a sketch of typical velocity 

distributions determined by this procedure. Plots of 2 ttV Z versus Z 

v v v 

were made assuming that the curve should asymptote to the value of the 
vortex strength once a distance had been traversed that was sufficient 
to account for all the vorticity. Thompson encountered some difficulty 
with his asymptotic method for calculating the vortex strength. The 
author found that the entire vorticity field had not been covered inside 
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the traversed area since the vortices were more diffuse than he 
expected. The author recognized this fact and stated that the 
given values for the vortex strengths were low for effective velocity 
ratios of 6 and 8. Thompson's method assumed the vortices were dis- 
crete and the velocity field fell off as 1/r, as in the single vortex 
filament. Since the vortices were diffuse, there was interaction 
between them which would also keep the asymptotic method from approach- 
ing the true value of the vortex strength. 

Fearn and Weston (ref. 1) conducted a study of the velocity field 
associated with a jet in a cross flow. The purpose of their study was 
to relate the velocity field, in cross sections perpendicular to the 
vortex curve, to the vortex properties through simple analytic models. 

A four inch diameter, unheated, circular jet directed normally through 
a 4 ft. by 9 ft. ground board was employed in the experiment. Velocity 
measurements were taken for effective velocity ratios from 3 to 10 and 
for a range of downstream distances of 2 to 45 jet diameters. The 
measurements were taken with a rake of seven yaw-pitch probes that was 
traversed in cross section planes. 

The authors presented two models for the contrarotating vortex 

pair associated with a jet in a cross flow. Like Thompson, the authors 

assumed that the vortex properties change gradually in the X direction. 

v 

In both models the vortices were treated as though they were two- 
dimensional and no attempt was made to account for an axial velocity 
component. The filament model approximated the vortex pair with two 
vortex filaments of strength +r and located at Y y = +h . The diffuse 
model assumed a Gaussian distribution of vorticity within each of the 



vortices. In both models the projection of measured velocities onto 
cross section planes were used to infer the location and strength of 
the vortices at that cross section. 

The filament model used measured velocities along the Z^ axis to 
determine the vortex properties. It was felt that the large upwash 
velocities (W^) , that occur along the axis would give the best re- 
sults in determining the vortex properties. The filament model assumed 
the velocity in a cross section was the result of the superposition of 
the free stream component of velocity in the plane of the cross section 
with the velocity induced by the vortex pair. Reference 1 gave the 
equation for the upwash velocity along the axis as 

W - rh/[?r(h^ + Z ^)] - U sin <f> (1) 

v v « v 

where $ is the angle between the Z and the Z^ axis. The vortex pro- 
perties were determined by fitting the measured upwash velocities in a 
least squares sense to equation (1). 

The diffuse model used a large number of measured upwash velocities 
in a cross section to determine the vortex properties. By fitting the 
measured upwash velocities in a least squares sense to the equation for 
the velocity predicted by the model, the strength , spacing h Q , and 
diffusivity 0 , of the two diffuse vortices were determined. Figure 5 
illustrates the coordinate system used in the development of the diffuse 
model. The vorticity u> was assumed to be 



8 


where m is the maximum vorticity of each vortex. The velocity at 
o 

any point in the cross section was assumed to be the result of the 
superposition of the velocity induced by the vortex distribution given 
in equation (2) and the component of the free stream velocity in the 
cross section plane. 


V = 


_o2_ 2 

(1-e 3 r l ) 


-S 2 r 2 
(1-e 15 r 2 ) 


2tt 


- U sin cj) e„ 
» v Z 


(3) 


The integrated strength , of a single diffuse vortex was defined to 
be 


T 

o 



2 it 

/ 


o 


r -B 2 r 2 , , 0 

a) e rdrdB 

J o 
o 


(4) 


Figure 6 shows the measured velocity vectors projected into a cross 
section plane together with the velocity predicted by the diffuse model. 
It is seen that the diffuse model provides an adequate description of 
these velocities. 

The authors were able to obtain several analytic relationships 
from the model. The effective strength V of each diffuse vortex was 
defined to be the net flux of vorticity across the half plane of the 
cross section and was given by 


tt/2 00 

T = / / o)(r,0) rdrdO (5) 

-n/2 o 


The effective spacing or center of vorticity h , was defined as 
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.. tt/2 00 

h = — J / o){r,6) Y rdrd0 (6) 

r j j v 

-tt/2 o 

Since diffusion and cancellation of vorticity were expected across the 
symmetry plane, T may be smaller then T . The authors pointed out that 
the effective strength and spacing determined by the diffuse model was 
assumed to be equal to the strength and spacing determined by the fila- 
ment model. Equations relating the parameters of the diffuse model were 
found by evaluating equations (5) and (6). The results were 


r = T erf (£h ) 
o o 


(7) 


and h = h /erf(3h ) 

o o 


where 


erf (3h ) = 
o 




( 8 ) 


is the error function. 

The relationship between the diffusivity and the vortex core size , 
was given by 


8 = 1.121/r c (9) 

where r was defined to be the distance from the vortex center to the 
c 

point where the maximum tangential velocity occurs. The authors presented 

most of the data in non-dimensionalized graphical form. The parameters 

h, h , and r were non-dimensionalized by the jet diameter. The 
o c 

strength of the vortices was non-dimensionalized by the quantity 2DU w 
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i.e. y - F/^DU^ and y q r^/2DU oo . The quantity 2DU^ was shown by 
Chang-Lu (ref. 7) to be equal to the roll up of the vorticity around 
a two-dimensional cylinder. 

The results of the experiment by Fearn and Weston (ref. 1) pro- 
vided several interesting implications. The authors found the vortex 
strength was essentially constant for each velocity ratio and could be 
described in a linear form 


Y 


o 


AR 


( 10 ) 


where R was the effective velocity ratio. By fitting the data for each 
velocity ratio separately, the constant was determined to be, A = 0.72. 
The fact that F q was related to R in a linear manner suggests that F 
was a function of the jet exit velocity and diameter. The authors 
presented a qualitative description of the vortex system: 

The vortex pair is formed very close to the 
jet orifice with an initial strength that is directly 
proportional to the speed of the jet at the orifice 
and to the diameter of the jet. The vortices are de- 
flected by the cross flow and they diffuse at a rate 
which is a function of the arc length along the vortex 
curve, but which is a weak function of effective velocity 
ratio. The vortices gradually weaken each other by 
diffusion of vorticity across the symmetry plane, (p. 1671) 


To summarize, Fearn and Weston were able to predict the velocity field 
induced by the vortex pair to an adequate degree by the diffuse model. 
Equations (6) through (9) were obtained relating several of the vortex 
parameters , 

It was later shown by Fearn (ref. 8) that 


- 1/2 


BD a (s/D) 


( 11 ) 
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where s was an arc length along the vortex curve and D was the jet 
diameter. This relationship was not without physical significance, 
in that equation (11) could be related to a "kinematic or eddy" 
viscosity for turbulent flow. 



EFFECT OF VARYING y IN THE DIFFUSE MODEL 

o 

The equation y = AR, given by Fearn and Weston, is an Important 
result of their diffuse model. Since the equation will be used in 
almost any effort to extend the diffuse model, a study is conducted 
as to its validity. An investigation is made into the procedure for 
determining the constant A, together with the effect that varying y 

o 

has on the vortex properties. 

The investigation is performed by varying the value of A in the 
two-parameter, diffuse model computer program developed by Robert Weston. 
The program is listed for reference in Appendix 1. The computer program 
utilizes the method of differential corrections (ref. 9), to fit the 
measured upwash velocities in a cross section to the Z component of 
equation (3). The computer program sets y^ equal to a constant and 
varies h^ and 3 for the best fit, in a least squares sense, to the 
measured upwash velocities. Once a fit is obtained for a cross section, 
the program calculates the ’’quality" of the fit in terms of a standard 
deviation o (ref. 10). In this study the computer program calculates 
the quality of the fit in terms of a percent standard deviation 


a = r- 5 — 100 

w V 

max 


( 12 ) 


where V x is the maximum upwash velocity induced by the vortex pair 

alone. Through an iterative procedure, the program is able to search 

and find the value of y = AR, which corresponds to the minimum a and 

° w 

thus the best fit to the upwash velocities. 
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In this study only the measured velocity data for an effective 

velocity ratio of 8 is used since it is the most extensively studied 

by Fearn and Weston. Figure 7 is a plot of o versus the constant A 

w 

in = AR. Table 1 contains information showing the size, location. 

and number of velocity measurements for each cross section in Figure 7. 

The overall curve is calculated by fitting the measured upwash 

velocities for all cross sections together. It is interesting to note 

from figure 7, that although a minimum point does exist for most of 

the cross sections, only a cutoff for small values of y^ is predicted 

by the diffuse model. The value of y^ can increase indefinitely 

(increase A) once a certain value of A has been reached, with no 

significant increase in a . To determine if this trend is a function 

w 

of the number of velocity measurements in a cross section, the measured 

velocities from Harms’ experiment (ref. 3) are input into the two- 

parameter, diffuse model program. Table 1 illustrates that one of the 

cross sections studied by Harms contains more velocity measurements 

than a similar cross section of Fearn and Weston. Figure 8 shows the o 

' w 

versus A curve that is the result of Harms’ measurements. It can be 

seen in figure 8 that the same large plateau exists in the a curve. 

w 

It appears that the number of velocity measurements does not change 

the large plateau in the a curve for increasing y 

w o 

In an effort to gain some insight into why the plateau occurs in 
the curves, the manner in which the diffuse model describes the upwash 
velocities for a range of A is investigated. The upwash velocities pre- 
dicted by the model are examined for values of A equal to 0.4, 0.72 and 
5.0. The cross section at X/D = 8.3 is examined since it is the most 



14 


extensively studied cross section for an effective velocity ratio of 
8. Figure 9 is a plot of the upwash velocities (W^) , versus Z^ along 
the symmetry plane of the cross section. The upwash velocities are 
non-dimensionalized by the free stream velocity. Figure 9 shows that 
the model describes the upwash velocities adequately for values of A 
equal to 0.72 and 5.0. It is also seen that the model cannot describe 
the upwash velocities for an A = 0.4. Figure 9 infers that a certain 
value of y is necessary to describe the upwash velocities and once this 
value of y is reached the net flux of vorticity across the half plane 
remains constant for increasing y q . This is also illustrated in figure 
10, which is a plot of y versus the constant A. For values of A > 1.0, 
it is seen that y remains essentially constant for each cross section. 

The plateau in the curves in figures 7 and 8 can then be ex- 
plained from 

Y = y (Bh ) (7) 

o.o 

together with the knowledge of how the vortex properties vary with y q . 

Figures 10, 11 and 12 illustrate that the properties y, h and r^ remain 

essentially constant for A. > 1.0. Figure 13 shows the vortex spacing 

h decreases rapidly with increasing y . The diffuse model attempts 
o o 

to keep the value of y in equation (7) constant by shifting the various 
vortex parameters. For increasing y^ , the term erf (fth^) in equation 
(7) must decrease to keep y constant. The properties r^_ or 3 remain 
essentially constant while h Q decreases in such a way that erf (&h Q ) 
decreases to keep y constant. For decreasing y Q , the term erf (gh^) 
must increase to keep y constant. The parameter r^ decreases 
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(3 increases), as the vortices tend toward filaments in an effort to 

match the upwash velocities. The vortex spacing h also increases, 

o ’ 

but the term erf (3h Q ) asymptotes to a value of one and cannot in- 
crease further. The result is that for y Q below a certain value, the 
effective vortex strength y decreases (figure 10). The vortices are 
then too weak to describe the measured upwash velocities which causes 

a to increase, 
w 

In summary, this investigation did not find an error in the value 
of A * 0.72 given by Fearn and Weston, since the overall fit for R = 8 
does have a physical minimum in this region (figure 7). This study 
does show that any value of A > 0.6 will work almost as well. 



■ EXTENSION OF DIFFUSE MODEL 


Fearn and Weston present most of the vortex properties in 
graphical form. It will be more convenient to an aircraft designer 
if the vortex properties are expressed algebraically as a function 
of arc length along the vortex curve for a given effective velocity 
ratio. This study, as stated previously, will consider only the 
diffuse vortex properties for an effective velocity ratio of 8. 

Fearn and Weston (ref. 1) present several equations relating the 
vortex properties, which are listed again for reference. 


y = y Q erf (Bh Q ) (7) 

h = h Q /erf (8h Q ) (8) 

Y q - AR (10) 

In addition, from reference 8 it is found, 

6D a (S/D) -1/2 (11) 


It is seen from equations (7), (8), (10) and (11), that 4 equations with 
5 unknowns ( y, y , h, h Q , 3 ), are available to describe the vortex 
properties. An additional equation is needed together with an explicit 
statement concerning equation (11) in order to describe the vortex 
properties algebraically. Before the projection of velocity field onto 
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a cross section plane can be reconstructed from the vortex parameters 
together with equation (3), a description is needed for the vortex 
curve. Ftarn and Weston (ref, 1) present an empirical equation for 
the vortex curve, which is given by 


Z/D 


b-y Ci 

a R (X/D) 
v 


(13) 


where a^ = 0.3473, b y = 1.127 and c^ = 0.4291. Values of <j> v and s/D 

can then be calculated for a range of X/D. 

In this study the two-parameter diffuse model program is used to 

generate the vortex properties at each cross section for an A = 0.72. 

From equation (9) , it is seen that a one to one correspondence occurs 

between 8 and r . Since the diffusivity 8 is more difficult to 
c 

visualize than the core size r , the latter will be given in the 

c 

figures. For convenience, 6 will be used in determining the vortex 
properties and converted to r^ through equation (9). Figure 14 is a 
plot of the vortex core size versus s/D. It is seen that no data is 
available near the jet orifice (s/D < 5), therefore, any empirical 
equations for r c or 8 will be an extrapolation in this region. This 
will also be true for any additional equations obtained from the re- 
sults of data presently available. From equation (11) » the simplest 
possible description of 8 is given by 


8D 


(s/D) 


1/2 


( 14 ) 


where B is a constant. 
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The constant B is determined by fitting the values of 3D predicted by 

the diffuse model, in a least squares sense to equation (14) . It is 

noted that equation (14) infers that the vortices approach filaments 

at the jet exit. Figure 14 shows that the empirical equation gives 

an adequate description of the vortex core size. Since 3 (or r ) is 

c 

a weak function of the effective velocity ratio, equation (14) will 
give a fairly reasonable description of 3 for other effective velocity f 
ratios as well. 

To make the system of equations complete, an additional equation 
is required which does not introduce still another unknown. One 
possibility is an equation describing the vortex spacing h , as a 
function of s/D. By examining the values of h^ predicted by the diffuse 
model, two possible options are obtained for the vortex spacing. Both 
options are based on assumptions for the starting positions of the 
vortices in the region near the jet orifice. Figure 15 shows the two 
options and 4 equations that this study will use to describe the vortex 
spacing. 

Figure 15 shows that option 1 assumes the vortices start as two 
concurrent vortex filaments at the jet exit. Since the contrarotating 
vortices are concurrent at the jet exit, this model forces the vortex 
strength to zero at the orifice. Two equations are shown that describe 
the vortex spacing adequately, 


or 


h 

o 

.D 


- s/D 
C (1-e ) 


h 

o 

D 


(15) 


s/D 

(s/D)+C2 


(15a) 
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where C, and are constants. 

Figure 15 shows that option 2 assumes the vortices start as 
vortex filaments that emerge from the side of the jet orifice. Option 
2 gives the possibility of a non-zero vortex strength at the jet ori- 
fice and will resemble the qualitative description of the vortices 
given by Fearn and Weston (ref. 1). The equations describing the 
vortex spacing are simple modifications of the two previous equations, 


(16) 
CL6a) 

The extended model consists of equations (7), (8), (10) and (14) 
together with one of the four equations for the vortex spacing. Each 
option of the extended model is denoted by the equation that is used to 
calculate h Q . The undetermined constants in the 4 options for the 
vortex spacing are determined by fitting the values of y predicted by 
the diffuse model, in a least squares sense, to the equation for y 
obtained from the extended model. Some question may occur over the 
decision to fit to y instead of the more straightforward method of 
fitting to . It is felt that the vortex strength is the most im- 
portant property and the best description of the strength is obtained 
by fitting to y . In addition, a later chapter will show that a large 
amount of uncertainty exists in the values of predicted by the 


h _s/D 

~ = C(l-e R ) + ,5 


or 


s/D 


D C 1 (s/D)+C 2 


+ .5 


where C, and are constants. 
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diffuse model. A description is given as to how the constant is 
determined in one of the options. For example option (la) will 
calculate the vortex strength y , from 


y = ARerf 


(s/D) 


1/2 


-s/D 
C(l-e *~) 


(17) 


The undetermined constant C is calculated by fitting the values of y 
predicted by the diffuse model, in a least squares sense, to equation 
(17). The constants C, and in the other options are found in 
a similar manner. 

Figure 16 is a summary of the extended models for an A - 0.72. 

The figure illustrates all the equations for each option together with 
values for all the constants. Figures 17 through 19 illustrate how 
well the different options describe the vortex properties for an A = 0.72. 
Figure 17 shows the vortex spacing , together with the empirical re- 
sults from the four options of the extended model. It is seen in 
figure 17 that there is little difference in options (la) and (lb) . 
Similarly there is little difference in options (2a) and (2b) . It is 
interesting to note, that for s/D greater than 10, it makes little 
difference whether the vortices start at the origin (options (la) and 
(lb)) or emerge from the side of the jet exit (options (2a) and (2b)). 
Figure 18 shows that the vortex strength y , is independent of initial 
vortex position for s/D greater than 10. It is also seen in figure 
18 that y does depend quite drastically on the initial vortex position 
for s/D less than 5. No conclusion can be reached on the value of y 
near the jet orifice until more data on the vortex strength and spacing 
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is obtained in this region. Figure 18 also shows the standard de- 
viation o of each option in fitting to y . Figure 19 shows that the 
effective vortex spacing h , is> independent of initial vortex position 
for s/D greater than 4. 

To illustrate how the model is used to calculate the vortex 
properties, a sample calculation will he made with the use of option 
(la). For a cross section located at X/D = 15.23, a value of <J> v = tan 
[d(Z/D)/d(X/D) ] can be calculated from equation (13). By a numerical 
integration of equation (13) , a corresponding value of s/D can be ob- 
tained. In this manner the values of 4> v = 18° and s/D - 20.5 are ob- 
tained. For an s/D of 20.5, the vortex properties are calculated as, 

Y q « (.72)8 = 5.76 

- 20.5 

h /D = 2 .04(l-e 8 ) = 1*88 

o 

6D = 2. 11/ (20. 5) 1/2 = 0.466 

y — y erf(8h ) - 4.537 
0 o 

^ ' 

h/D = (h Q /D)/erf (6h Q ) = 2.396. 

The projection of the velocity field onto this cross section plane can 
then be calculated at any point throughout the cross section from 
equation (3) . 

In summary, several options are obtained to extend the diffuse 
model given by Fearn and Weston. With the use of these options, the 
vortex properties and the projection of the velocity field onto any 
cross section plane can be calculated. It is noted however, the 
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options are extrapolations of available data in the region of the jet 
orifice. For s/D > 10, the vortex properties are adequately described 
by any one of the four options. 



UNCERTAINTY IN EXTENDED MODEL 


Since the extended models are the result of curve fitting to the 
diffuse vortex properties, some insight into the uncertainty of the 
diffuse vortex properties is necessary. The criteria that are used 
to set the limits on the uncertainty of the diffuse vortex properties 
are established with the use of figure 7. As stated previously, the 
lower cutoff on y q is distinct. For each cross section in figure 7 
the lower or percentage cutoff is defined to be the value of the con- 
stant A that corresponds to twice the value of a w at t ^ ie ^itLf mui11 
the cross section. The plateau on the c w curves raises some question 
as to a proper cutoff limit for large values of A. Since the constant 
A can increase indefinitely without a large increase in the value of 

a the value of A = 2.0 is arbitrarily chosen. This value of A is 
w 

sufficiently far out on the plateau that the vortex properties (with 
the exception of h Q ) are not changing rapidly with A. With this 
criteria established, figures 10 through 13 are used to determine the 
diffuse vortex properties at the percentage (lower) and plateau (upper) 
cutoffs. It should be noted that the percentage cutoff will vary from 
cross section to cross section but the plateau cutoff remains at 
A - 2.0. 

Figure 20 is a plot of the core size from the diffuse model for a 
constant A = 0.72. The uncertainty bars are shown for each cross 
section. The bars marked with the double tick mark are the result of 
the percentage cutoff for each cross section. The bars marked with a 


23 



24 


single tick mark are the result of the plateau cutoff. It is 
necessary to determine if the empirical equation for r^, equation (14) 
together with equation (9) , will still describe the vortex properties 
throughout the region of uncertainty. Equation (14) is fit in least 
squares sense to the values of 3^ predicted by the diffuse model for 
constants of A equal to 0.48 and 2.0. These two values of the constant 
A correspond to the percentage and plateau cutoff as determined from 
the overall fit for all cross sections in figure 7. These two values 
of A will then be used to set uncertainty limits on the extended model. 
Figure 20 shows that the empirical equations will describe the vortex 
core size adequately in the region of uncertainty. 

In the remaining portion of this study, options (la) and (2a) are 
used because of their relative simplicity. The diffuse vortex pro- 
perties shown in figures 21 through 26 are determined from the diffuse 
model for a constant A = 0.72. Figures 21 through 23 show the uncer- 
tainty in the diffuse vortex properties together with the empirical 
descriptions given by option (la). Figure 21 Illustrates that there is 
a large uncertainty in the vortex spacing h Q given by the diffuse model. 

It is also seen in figure 21 that option (la) still gives an adequate 
description of the vortex spacing in the region of uncertainty. Figure 
22 shows the uncertainty in the effective vortex strength y , from the 
diffuse model together with the filament model results (ref. 1). 

Extended model (la) is seen to describe the vortex strength adequately 
and the standard deviation a, in fitting to y is given for each curve. 
Figure 23 shows that there is very little uncertainty in the effective 
vortex spacing given by the diffuse model. For all but one cross section, 



25 


the uncertainty bars remained within the symbol width. There is a 
discrepancy in that some of the filament model results lie outside of 
the uncertainty limits of the diffuse vortex properties and the extended 
model curves. It is believed that at these points an insufficient 
number of velocity measurements were taken for the filament model to 
adequately describe the upwash velocity distribution. This will cause 
the filament model to give unreliable values for the vortex properties. 

Figures 2A through 26 are similar to figures 21 through 23 except 
the empirical curves are determined by option (2a). It is seen in 
figure 25 that option (2a) encounters difficulty in describing the 
vortex strength y, for a constant of A = 2.0. 

In summary, the extended model will give adequate descriptions 
throughout the range of uncertainty in the diffuse vortex properties. 

The one exception is that option (2a) cannot describe the vortex 
strength for a value of A = 2.0. 



SUMMARY AMD CONCLUSIONS 

An investigation is conducted into the diffuse vortex model 
given by Fearn and Weston (ref. 1). The equation, y^ - AR, presented 
by the authors, is examined in detail as to the procedure for determin- 
ing the constant A. As in reference 1, the value of A = 0.72 is found 
to give the best description of the upwash velocities in a cross 
section for an effective velocity ratio of 8. However, this study 
has also shown that, in a practical sense, any value of A greater 
than 0.6 will work almost as well. 

The diffuse vortex properties presented by Fearn and Weston, are 
extended from a graphical to an algebraic description of the vortex 
properties for an effective velocity ratio of 8. The extended model 
consists of the analytic equations given in reference 1, together with 
empirical equations given in this paper. With the use of the extended 
model, the vortex properties together with the projection of the 
velocity field in a cross section plane can be calculated for any 
cross section in the flow. However, in the region near the jet orifice, 
the extended model represents an extrapolation of available data. The 
extended model gives an adequate description of the vortex properties 
for s/D > 10. The uncertainty in the results of the extended model is 
also investigated. 
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Figure 3. Sketch of Jet Centerline and Vortex 
Curve with Coordinate Systems 
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Figure 8- Standard Deviation versus y , R=*8 
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Figure 13. Variation of Vortex Spacing 
with , R « 8 
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EXTENDED VORTEX MODEL 

Y = AR where A = 0.72 
o 

3D = B/(s/D)~ 1/2 where B - 2.11 (r = 1.121/3) 




Y = Y 0 erf (6h o ) 

h = h /erf (3h ) 
o o 


Figure 16. Summary of extended model, A ~ 0.72 
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Figure 21. Uncertainty in Vortex Spacing, 
Option (la) , R = 8 

















TABLE I 

SUMMARY OF CROSS SECTION DATA 


VELOCITY 

X/D 

Z/D 

s/D 

<f>v 

deg. 

CROSS SECTION SIZE 

MEASUREMENTS 




JET DIA. 

EFF. VORTEX SPACING 

49 

2.13 

4.87 

5.52 

51.7 

2.50D x 2.48D 

2. 17h x 2. 16h 

135 

5.24 

7.56 

9.65 

33.0 

4.99D x 6.97D 

2. 72h x 3. 80h 

36 

6.01 

7.98 

10.53 

30.7 

2.49D x 2.50D 

1, 48h x 1.49h 

216 

8.34 

9.14 

13.11 

26.2 

6.28D x 8.44D 

3.02h x 4.05h 

140 

15.23 

11.83 

20.5 

18.0 

4.49D x 6.52D 

1. 86h x 2 . 70h 

49 

35.57 

15.85 

41.29 

10.9 

2.94D x 2.99D 

0 . 87h x 0 . 89h 

t>340 

15.0 

11.16 

20.6 

17.7 

5.41D x 8.43D 

2 . 36h x 3.68h 


Fearn and Weston (ref. 1) 
Harms (ref. 3) 












APPENDIX I 


TWO-PARAMETER DIFFUSE MODEL COMPUTER PROGRAM 




c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 

c 

c 


c 


THIS PROGRAM FITS A GAUSSIAN- DISTRIBUTED VORTEX PAIR TO THE CROSS-SECT I CN 
VELOCITY FIELD FOR A JET IN A CROSSFLOW 
ALL VELOCITIES ARE IN FT. /SEC. 

A DENOTES THE CIRCULATION IN FT**2/SEC ANC GAMMA DENOTES THE DIMENSIONLESS 
C I RCUL AT ION, GAMM A= A./ ( 2*C*V INF ) , HHERE D IS THE DIAMETER OF THE JET AND 
VINF IS THE FREE STREAM VELOCITY. 

0 DENOTES THE DIFFUSION CONSTANT IN THE GAUSSIAN DISTRIBUTION. IT HAS UNITS 
OF 1/FT. 

C DENOTES THE HALF SPACING BETWEEN THE CENTERS OF THE TWO DIFFUSE VGRTICESi 
H DENOTES THE HALF SPACING IN JET DIAMETERS, H=C/C. LET BETA DENOTE THE 
DIMENSIONLESS DIFFUSION CONSTANT, BETA=B*C. 

DIMENSION RMSGt 6 ) ,RKST< 6 ),P (6 ) , St 3 ) ,MBC (5) 

COMMON I YZCOD, I Y Z , NR , NC ,A*B,C,ZC,VINF»VSINA,VYt2l,I2)*VZl2l,l2I, 
lY(3,I2),Z(21),RMSfNCRPF, I PR NT , NPTS, IRMS , RMSVZOt 2 , 5 ) , IPRTMiNRAK 
CCMMCN/ONE/MK,NCATA, IRUN t A) , R , AM JET , XA , Z A , SO t ANG VRT , VRANG 
RE A L*8 P,S,0 


2 


1 

IF 

IF 

IF 


IC 


DC 2 1=1,6 
RMSTt I)=0.0 


0TR=1.7A5329E-2 

READ 1, CONST, GPER, IPRTM, IRMS 

FORMAT <1 AX , F8.3, 17X,F8.fc,7X, H,6X,11> 

I PRT M = Q , NO MATRICES PRINT GUT. 

IPRTM-2 , VORT IC I TY MATRIX CCESNT PRINT 
I R M S = C * OVERALL FIT IS IN PER CENT 
GPER=-01*GPER 


IF I PRTMs 1 , ALL MATRICES PRINT GUT. 
OUT, BUT OTHERS DO. 

IF 1 R M S = 1 OVERALL FIT IS IN FT/SEC 


CC 10 1=1,6 
PMSG t I ) = C*0 
NG = C 


N B C = C 


20 N = 0 

V INF = C • 

AM JET=0. 

CONSULT FORMAT NC. 50 FOR CESCRIPTICN OF VARIABLES IN FOLLOWING STATEMENT 
REAC(5,30,ENC=9S59) ANG VRT , ZO ,NR AK , NDRPF , NDRPL ,NC ,G AMMA , BET A,H,NZ, 
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1ZOC 

C NZ IS THE NUMBER OF THE RAKE CENT A INING THE VORTEX CENTER. 

C ZOC=STEF CHANGE IN 1 USED IN ZC PERTURBATION ANALYSIS 

30 FORMAT! lOXt F4.1, 8X, F 5.2, 5X, 11 , l 2X, IX , X7X , 11/ 12, 37X,3 IF5.2* 6X1 f I 1* 

I F 7. 2 > 

C AX ANC Cl USED IN REFERENCING RAKE POSITIONS IN THE CROSS- SECT ION 
40 AX- T AM AN GVRT*DTR > 

C 1= SORT ( 1. + A1*A1 ) 

C NR=NC . CF ROWS OF CATA 

NR = 7 *NRAK— NDRPF-NDRPL 
11-7 

IF ( NZ- EQ.NRAK) IZ=7-NCRPL 

WR I T 6 ( 6* 50) ANGVRTf ZOr IZ ,NZ , ZOC , NRAK ,NCRPF , NORPL *NC , GAMMA $8ETA ,H 
50 FORM AT ( 1 1 INPUT CATAC VORTEX ANGLE sSF5.lt* DEGS7X, ' ZO *SF6.2f 
1* IN. FROM PR0BESl2t # OF RAKE* « 12 » 7X t ’DELTA ZO *SF5.2/ 

X 1 5X» lit* RAKES OF INPUT S12»* PROBES FROM FIRST RAKE ARE D 

2EIETED*/33X, IX, * PROBES FROM LAST RAKE ARE DELETED* /15X, 12* * YAW 
3PC$ I T IONS * // * INITIAL ESTIMATESC G APR A =*,F5.2/ 21X, 'BETA = S 
4F5.2/ 21 X* 1 H = S F 5 . 2 ) 

R = G * 

C DC LOOP 110 READS IN VELOCITY FIELD DATA ANC RAKE POSITIONS* AND CONVERTS 

C THESE TC THE CROSS-SECTION COORDINATE SYSTEM 
DC 11C 1=1, NRAK 

READ (5,60 I RUN ( I ) , ANGR AK , Q ,RHO , X4, Z4 , A 2 , ( Y (I , J ) # J= X *NC ) 

60 FORMAT (3X, 13 »11X,F4. l>3X,F6.3,5X,F8.6t2l4X,F5.2) , 3X ,F5. 2/12F6. 2 ) 

R = R + A 2 

IF ( I. EC. 1 ) Cls*Z4-Al*X4 
04= { Al*X4— Z4*C 1 ) /C X 
VINF=VINF+SQRT( 2.+C/RH0 ) 

CA~lANGRAK-ANGVRT>*OTR 
$ INC A = S I MDA ) 

COSC A = CO ST DA ) 

L I - I 
. L=N*1 
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JR QUALITYi 


M= L 

IF ( I . EC. 1 ) Ll^UNDRPF 
IF < I .EG* 1 ) L~L1 
N = N + 7 

IF< UEC.NRAK) N = N-NCRPL 
I F ( I • N E • N Z) GG TG 7C 
XC=X4/4. 

ZD-Z4/4 • 

R AKEZC^ZC 
VRCC S=CO$DA 
VRAl\G=ANGRAK*DTR 
I Z-N 

Z0= ZO*COSDA 
ZOC=ZCC*CGSOA/12. 

70 CG 9C K= 1 »NC 

REA 0(5,80 <VYCJ«K) 9 VZ{ J«K) »Z(J ) I 

80 FORMAT! 11F7*2/1CF7*2) 

IF ( I * EG* NRAK .AND. NCRPL*GT.3) READ(5#80J DA 

CG 9C J=UN 

VY I J~L 1 + 1 > K ) =VY ( J * K ) 

90 VZ< J-Ll+ltK)=Z(J)*SINCA + VZ ( J,K )*CCSDA 
DC ICC J=L,N 

100 Z< J-Ll+1 )=044< J-L+X U*2.*CCSCA 
110 IF(UEe-l) N^N-NCRPF 
R=R/NRAK 
GAMMA=CGNST*R 
I F ( N Z» EQ * i ) IZ=IZ-NCRPF 
V I NF -V INF/NR AK 

C COMPUTE NOMINAL JET MACH NUMBER USING ESTIMATE CF PINF 12130; PSFA) 
AMJET=R*5CRT (C/(C#7*213C* I ) 

NCMR- 1 NT ( R+0 * 5 ) 

CALL CMAJINGMR* A 1* A2, A3 ) 

NK = 0 

IF! A 8 $ ( G ** A 3 ) *LT .2* ) MK-3 



Qg 
m g 



IF(^ES(Q-A2)*LT*2* ) MK=2 
IF ( AeS<Q-Al).LT.2.) MK=1 
VSINA=VINF*S IN! ANGVRT*DTR ) 

CG 12C J = It NR 
120 Z<J)«Z(NR)-Z(J) 

20- Z ( I Zl + Z0 
R EFZO-20 

C DETERMINE INITIAL ESTIMATES 
A=2.*VIN F*G AMM A / 3 • 

C=H/3. 

E=BETA/C 
BS A V = B 
CS A V=C 

*RITE(6, 130) Af£iC»VINF*R* ( IRUN (I ) , 1= 1 t NR AK } 

130 FORMAT ( / 14H A { FT**2/ SEC ) = * F 1C . 2 , * B t 1/FT ) = 1 » F 10 .4 * ' G < FT I * • *F10 .3/ 

1/' AVG VINF ^ , ,F6*1/* AVG ft =*,F6.2/' RLNS READ ING*,4f4) 

PRINT 135 * AM JET » MK 

135 FORM AT ( / 1 NOMINAL MACH NO. - * t F5 .3 f 1CX » • G CONDITION = ' ♦ 1 2 ) 

A1=VRANG/DTR 

PRINT 14CtNZ*XC,ZC, A1 

140 FORMAT!/' PROBE 4 OF RAKE', 12,' AT X/C =',F6.2,' f Z/0 = HF6.2|2X* 

1* RAKE ANGLE = ' * F 5 * 1 » 1 DEG * ) 

CALL AUXIL $$$$$$** 

C SET Z SPACING OF 5 VORTEX CENTER CASES. FCR PERTURBATION ANALYSIS. 

RMSVZCf 2, 1)=0. 

RMS V ZC ! 2 * 2 ) =-24 . * ZQC 
RMSVZO(2,3)=-ZOC*12. 

RMSVZC(2i4)=Z0C*12. 

RMS VZC ( 2 » 5) =24* *ZCC 

ZOSAVZ=ZO 

RMS=0* 

MG=- 1 

C MG IS THE [ NOEX OF THE GAMMAO AND VELOCITY COMPONENT BEING CONDUCTED. 

I GF = C 


Ui 

CO 




CO 720 KG = 1 » 3 

C DO ANALYSIS FOR 3 SETTINGS CF G A MM AO 
IFIKG.EQ.2) AM 1 • FGPER ) * A 
IF ( KG • EG • 3 ) AM l.-GPER)*A/( 1 * FGPER ) 

ZOS A V = ZOSAVZ 

C DC ANALYSIS FOR FIT TO Z VELOCITY COMPONENT. 

JLZ 35 1 

C IF LAST TIME BOMB EC CUT, RESET INITIAL GUESSES* OTHERWISE, START 
C RESULTS CF PREVIOUS TIME. 

IFIRMS.EG.O. > GC TO 34C 
BSAVZ-B 
CS AVZ-C 
GC TC 350 
340 BSAVZ= BS AV 
CS A V Z-CS A V 
350 ZO=ZCSAV 
I YZ -3~ JL Z 
IYZCCDMYZ 

C IYZCGCEM FCR VY CALCULAT ICN -2 FOR V Z CALCULATION 
NIT s O 

C NIT IS NUMBER OF RE-SHIFTS CF ZC. LIMITED TC 9. 

IF ( KG • EG .1 ) .WRITE<6» 360 ) 

360 FORM- AT 1 1 HI ) 

370 MT«MTM 
Z0P“ZC*12. 

IFIKG.EQ.T) PRINT 38Q,NIT,Z0P 
380 FORMAT ( / ' ITER ZC =*,Ffc.2I 

IPRNT-0 
IR0 = C 

IF ( N IT *£G .9 ) GG TO 540 
C DC tCCP A 30 PERFORMS FIT FCR 5 ZO CASES 
CC 4 30 LZ-1 1 5 
I F ( L 2 . EQ *1 > GC TC 390 
IFILZ.EQ.2) Z0 5S ZG-2.*ZCC 
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C 

C 


c 


IF(LZ.EQ.3I ZO=ZC+ZCC 
IF t L Z • EQ -^« ) Z0=ZC+2.*ZCC 
IF ( LZ • EQ«5 ) ZO*ZC+ZGC 
390 [FIJLZ.EC.l • A N C * KG.NE.l) GC TO AOO 
IF{RMS.EG*0. 1 GC TO AlC 
esAv=e 
CS A V = C 
GC TC A2C 

ago 8SAv=eSAvz 

CS A V*CSAV Z 


aio e^esAv 
C-CSAV 

A20 CALL DIFCCR 
CIFCCR CALCULATES 


LEAST SQUARE FIT TC THE VELOCITIES 


RMSVZOU »LZ)=RMS 

IF ( RNS.EG.O. ) IRC-IRO^l 

IRC IS THE NUMBER CF BOMB-CUTS CF THE 5 ATTEMPTS • 


A 30 CCNTINUE 

Z0 = ZC-2* *ZOC 
I F t IR0.GT.2) GC TC A70 

IF(KC.EQ.l) PRINT AAO» t { RMS VZO ( 1 ♦ J ) ♦ I - X »2 ) * J = l * 5 ) 
AAO FOR N AT < 1 RMS ^'fFS^A, 1 AT VZC *»*F6.2} 

FI I ) s 0 *0 
P(2)=C.O 
P { 3 ) =0 • 0 
P (A ) “C .0 
P ( 5 ) *0* 0 
S( 1 ) =0 *0 
SC2)=C,0 

S(3)=C.O • . 

PERFCRM SUM MAT I CNS FCR PAR A 60L 1C FIT TC RMS VERSUS ZO 

CC A60 1*1,5 

IF ( R K S V Z C ( 1, I) • E C • C • I GC TC AfcO 
A 1 = 1 • 
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C 

c 


CC 450 J*l,5 


450 

460 


470 

480 


IF 


P( J)=PIJ)+AL 
IFtJ.GT.3I GC TC 45C 
SU) = S(J>+Al*RMSVZOtUl ) 
Al-Al*RMSVZ0t2f I > 


COM INUE 

CALL MINV3C P»C) 

IFtC.NE.G.) GC TC 560 

PRINT 480 » ( ( RMSVZOl I»J)»I = i»2i*J = lt5) 

FORMAT I / * NO PARABOLIC SOLUTION TC RMS VERSUS ZC *//’ 
U' R MS = • t F8 *4 * • AT VZO ='*F6.2>) 

IF(NIT.EG.9> GC TO 53C 
IF ( IRC.NE.5) GC TC 51C 


INPUT W A S 0 • / 


ALL 5 ATTEMPTS FAILED* LCCP 500 SEARCHES FOR A SOLUTION 8Y VARYING ZC UP 


AND DOWN. 

CC 5CC J=l,6 
IFU.EQ.l) ZC*ZC-3.*ZCC 
IF( J.EG.2) ZO=ZC-ZOC 
IFU.6Q.3) ZC=ZC-ZOC 
IFtJ.EQ.AI ZC-ZC38.*ZCC 
IFU.GE.5l ZC=5ZC + Z0C 
E~6$AV 


C -C S A V 
CALL CIFCGR 
ZOP=ZO*12. 

PRINT 49C,RMS»ZCP 

49C FORMAT!* RMS= , »F8.4,' AT ZG='tF6.2) 

IF(RMS.NE.O. > GC TO 37C 
5C0 CONTINUE 
GO T C 5 3 C 

C520 AN ATTEMPT THAT WORKED IS USED AS THE NEW CENTER CF THE ZO SPREAD. 

510 CC 52C J = 2 » 5 

IF(RMSVZOtl» JI.EG.O. ) GC TC 520 
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IF(J*EG*2> ZG=ZC-2.*ZCC 
IF(J.EQ.3) ZC=ZC-ZOC 
IF ( J *EQ. 4 > ZO=ZC*ZOC 
IFU-EQ.5) ZC=ZC*2.*ZCC 
GC TC 37C 
520 CONTINUE 
530 MG=MG*2 

f F ( JLZ *EC* 1 > IGF= l 
GG TC 72C 
540 PRINT 55C 

550 FORMAT t / * ITERATION LIMIT EXCEEDED IN FITTING RMS VERSUS Z0*1 
GC TC 47 C 

C COMPUTE NEW ZO (AT MINIMUM CF THE LSQ PARABOLIC FIT) 

56C CZ0=(S<l>*P(2)+S(2)*P(4HS(3)*Pl5)>/t24.*(S<l)*P{3)+S<2)*P<5H- 
1S<3)4P(6>)) 

A 1 = 5 . *ZOC 

IF( AES(DZO) .GT.A1) CZC= S I GN ( A 1 ? CZ C ) 

Z0= ZO-DZ C 

C METHOD IS ITERATED UNTIL THE ZC VERSUS RMS CURVE HAS ITS MINIMUM WITHIN THE 
C ZO RANGE. 

IF(AfiS(DZO).GT. ( 1 * 2*Z OC ) ) GC TO 370 
Z0P=Z0*12. 

PRINT 57C »ZOP 

570 FOR M AT l / ' FINAL VALUE FCR ZC **,F7.2) 

Al = „25*(RAKEZO+(ZOP-REFZC)/VRCOS-6. ) 

X4=XC~A1*$IN(VRANG) 

Z4=ZC+Ai*CCSCVRANG) 

Al = Z4-.347308*R**l. 126536*X4**.429137 
SC = ARCLNG(X4*R) * A 1*S IM ANG VRT*CTR ) 

PRINT 580fX4tZ4tSC 

580 FORMAT ( / 1 VORTEX CENTER AT X/C =*,F6.2t*t Z/C = * »F 6 • 2# 1 CK ♦ * ARC LEN 
1GTH (S/0) = 1 f F6 *2 ) 

IPKN T =1 

C NCW PERFORM LEAST SCUARES FIT AT THE MINIMUM RMS LOCATION. 
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590 CALL DIFCOR 

IF I KG .NE • 1 ) GO IC 71C 
I F { IPRTM.EQ,C) GC TC 71C 
CALL PRTMAT 


$$**$$$$ 


710 NG=MG+2 ' „ ^ „ 

C PERFORM SUMMATIONS FOR OVERALL RMS COMPUTATIONS FOR ALL CASES AT THE SINGLE 

C VELOCITY RATIO BEING CONSIDERED* 

PMST (MG) =RMS*RMS*NDATA 
I F ( J LZ * EC *2 ) GC TG 72C 
ZGSAVZ-ZC 

IFIRMS.EC.O. ) ZCSAVZ-ZCSAV 
720 CONTINUE 

IF I I C*F *NE * 1 > GG TC 730 

C IGF=1 IF ONE GF THE GAMMAC CASES BOMBED OUT IN FITTING Z COMPONENT • 

C THIS RUN IS THEN NOT INCLUDED IN THE OVERALL RMS CALCULATION* 

NBC=NEC+1 
M BO I NEC ) - IRUM 1 ) 

GC TC 20 

730 DC 7 AC 1=1,6 
7 AO RMSG ( I) = RMSGU ) 4RMSTT I) 

NG=NG+ND AT A 
GC TC 20 

5559 CC ICCGO 1=1,6 
1C0C0 FMSG ( I ) = SGRT I RMSG ( I ) / N G ) 

A 1= { U+GPER ) *CCNST 
A2= ( 1*-GPER ) *CGN ST 

PRINT 10C10, RMSG ( l ) ,RKSG l 2) , CONST , RMSG ( 3 ) ,RMSG ( A ) , A1 , RMSG ( 5 ) , 

1RMSG ( 6 ) , A2 

1C 010 F C R H A T ( / / / 1 X , 2 8 < 1H * ) / 2 9 F *RMS VALUES FOR OVERALL F I T*/lXfc 28 ( 1H* ) // 

I* Z RMS =*,F8*A, 1 , Y RMS-. s, »F8*At* FOR CENTRAL CCNST «#,F6*3/ 

2* Z RMS = * , F8* A, 1 , V RMS = * »F8. A , 1 FOR HIGH CONST = * »F6* 3/ 

3* Z RMS = * , F8-A , • , Y RMS = , ,F8,A, 1 FOR LCV CONST s’ 9 F6«3) 

A i=GP£ R*CCNST 

A 2 = (RMSG(3)-RMSG<5) )/(2 .*A1) 


cn 

u> 




A3 = (RMSG(3)+RMSG(5)-2.*RM5G( 1) )/<2.*Al*AlI 

A1=-C.5*A2/A3 

e=CCKST+Al 

A=RMSGI l >*A1*I A2 4A3*A1 ) 

IF ( IRMS. EC.Q ) PRINT lC02C,A f B,NG 

IG020 FORMAT I / * MIN IMUM RMS = # ,F8.<i, # PER CENT AT CONST = f >F6 . 3 // 1 5% 1 PC 
II NT'S CF CAT A* ) 

IF! l R M S • E G . 1 ) PRINT 10021, A, B,NG 

IC02I FORMAT!/* MINIMUM RMS =*,F8.4 t * FT/SEC AT CONST = * , F6 * 3/iT 5 » 4 PC I N 
ITS CF CATA* ) 

IF(NBC.EC.O) STCP 

PRINT 10030, (MBCU), 1*1. NBC) 

ICO 30 FORMAT!/* RUNS NCT USED'C 1 *5 ( 14, • , *)) 

STCP 

ENC 

SUBROUTINE CIFCCR 

C THIS SUBROUTINE PERFORMS LEAST SQUARE FITTING USING THE El FFERENTI AL 
C CCRRECTICN METHCC FOR A GAUSSIAN -DISTRIBUTED VORTEX PAIR, 

C ThC PARAMETER FITC GAMMA IS HELD CONSTANT e AND C ARE VARIED FOR BEST FIT. 
COMMON I YZCGC*I YZ ,NR »NC , A ,8 *C ,ZQ ,VINF , VS INA » VY 121,12) ,VZ(21,12) , 
lYl3,12>,Z12i),RM$,NCRPF, IPRNT,NPTS, IRMS tRMSVZO I 2 , 5 ) , I PRTM *NR AK 
COMMCN/ONE/MK,NCATA, I RUM 4 > , R , AM J ET » X4 , Z4,SC, ANGVRT , VRANG 
CCMMCN/DVC/F IB, F IC 
CATA ERRCR/4.0E-4/ 

C ERROR = ACCEPTABLE ERROR FOR VARIABLE COEFFICIENTS. 

IF! IPRNT.EC.2) CC TC ICC 

ITEP=C 

IYZCCCMYZ 

C TWC VARIA8LE COEFFICIENTS I B,C ) 

C 

C PERFORM SUMMATIONS FOR COEFFICIENTS IN NORMAL EQUATIONS 
10 SF82=0.0 
SFC2 = C * 0 
SFBC=0*0 


o> 
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c 


c 


SFER=C .0 
SFCR=G.O 
SR 1 2=0 *0 
NPT S = C 

CC 2C 1 = 1 , NR 
K=( 14NCRPF+6J/7 
CC 2C J= 1-tNC 

IF ( VY < If J") * GT • 5CC # ) GC TC 2C 
IF ( IYZ.EC.n V Z Y = V Y ( I»J) 

IF (IYZ.EG.2 ) V Z Y = V Z ( I * J HUS IN A 
RI=FIY(K»J)*Z(I) )— VZY 
F I A = CF ( Y l K» J > * Z ( I ) ) 

SRI2=(SRI2+RI*RI) 

SF62=SFB2 + F I B * F I B 
SFC2=SFC2 + FIC*FIC 
SFBC=SFBC + FIB*FIC 
SF ER = SFBR - FIB*RI 
SFCR = S FCR - FIC*Rl 
NPTS=NPTS+i 

20 CONTINUE 

RMS = S CRT { SR 12 /N PTS I 

SOLVE NORMAL EQUATIONS FOR THE CORRECTIONS TC THE ESTIMATES 
C =SFB2*5FC2 - SFBC*SFBC 
IFIC.EQ.C.) GC TC SSS 
E 1= l S F BR *SFC2 - S FCR*SF EC I /C 
C 1= ( $ F82* SFCR - SFBC*SFER)/C 
I F ( ITER.GT.2CI GC TC AC 
I F ( ITER. GT. 1C) GC TO 7C 
IFIAESIB1/B) *jGT.. 05) GO TO 30 
IF(AESICi/C) *LT . «C5 ) GC TC TC 

DAMP THE CORRECTIONS — : 

30 I F ( ITER.GT.5) GC TO 5C 
FRAC={ ITER+3I/2C. 

GC TO 60 


CTv 
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40 FRAOO.5 
GC TC 60 

50 F R A C = ( ITER—1 1/10 • 

60 B1»FRAC*B1 
C1=FRAC*C1 

C UPDATE ESTIMATES AND TEST FCR CONVERGENCE 
70 ITER= ITER + 1 

E=e+ei 

C=C+C1 

IF ( I PRNT .EG* 3 > PRINT 80 * I TER ,RM S , A* 8 »C 
80 F0RMAT(2GX* 1 I7ER = f f 13# 1 * RM S = • , F8*4 , • , A=«,F8.4,', 8^VF8i4,*% C** 
1 * F8 • 4 ) 

IF ( ITER. GT. 30 GC TC S9S 
IFI8.LT. C.) 0 =■— 8 

IF(A8S(B1/B).GT.EFRCR. CP. ABStCl/O.GT. ERROR) GO TC 10 
SR I 2 = 0 • 

C CCMPUT E STANDARD DEVIATIONS. 

CO 9C 1 = 1, NR 
K=i ( I4NCRPF+6J/7 
DC 9C J= l ,NC 

IF(VY(I,J).G7.5CC.) GC TC 9C 
I F { I Y Z • E C* 1 ) VZY*VY(I»J> 

IFUYZ.EC.2) VZY = VZ ( I,J)+VSINA 
RI = F (Y(K,J> ,Z( I) V Z Y 
SRI2=ISRI2+RI*RI ) 

90 CONTINUE 

R KS = SGRT (SRI 2/NPTS ) 

I F 1 IPMS.EQ. 1) GC TC 55 
I F ( I PRNT • EQ * 1 ) GC TC ICC 
C CALCULATE PERCENT RMS ERROR 

RM$=314# 15926*RMS*C/ ( A* (l.-EXP l -1 8*C )**2 ) ) ) 

95 I F { I PRNT * NE • 1 ) RETURN 
ICO GAMMA 0=1 .5*A/VINF 
B8AR=E*C 


crv 

o 






FD = 3 **C 

£RFe=ERF<EBAR> 

GEFF=GAMMAO*£RFE 

FtEFF^HD/ERFB e ^ fIlll 

VMAXUP=2.*VINF*GAMMAC*( 1 *-EXP (-BBAR**2 ) )./ CH0*3. 141593 J-V5INA 

IFC IPRNT.EQ.2) GO TO 14C 
E1*ICC.*RMS/ (VMAXUP + VSINA ) 

C PUNCH 

A1=8BAR/HC 
ei=,ci*Ei 

C1=VPANG*57. 29576 „ , 

PUNCH 10 5 . MK »NC AT A » IRUN ( 1) . R • AM JET , X4, Z4*hEFF *GEFF , 6 L*SD%ANGVRT ,CH 

1»VINF*GAMMA0*HD»A1 

105 FORMAT U 1*1 HO » 13* 14 , F 5 . 2 . F5. 3 1 2F5. 2 1 2F5* 3 ,F5* 4, 2F5* 2 »F4. 1 i F5*l 
l,F5.2,F5.3 f F6.4> 

C 4 $ $4S$ $ $ 

WRITE (6, 110) ITER t Al,ei»Cl,A f 0,CtRMS»El 
HO FORMAT l / / 1 FINAL RESULTSC I TER* ' » 1 3 . * » A1.B1.C1 -* *3F9*6// 

1 18X » • A= • . F9* 4 * * » B=»»F9.4 F ' f C- % F9 .4// i6Xt • R^S CF FIT =* F9«4 t 
2 * FT/SEC -».F6*2.9H PER CENT) 

IFt IRMS.EC.O) R M S =; E 1 

C SWITCH Y ZCCDE FCR ALTERNATE RMS CALCULATION 
IYZCCC=3-lYZ 
SR 12-0 • 

CC 12C 1*1. NR 
K=( I4NCRPF + 6J/7 
CC 120 J= 1 . NC 

I F < V Y < I.J).GT*5C0*) 00 TC 12C 
IF I IYZ*£G*2) VZY=VY( WJ i 
IF ( I Y2« EC * l ) VZY-VZt I.J MVS INA 
PI = F (YtK, J) * 2 ( I ) i~VZY 
SR 1 2 = ( SR 1 2 + R I *R I) 

120 CONTINUE 

SRI2=SCRT (SRI2/NPTS) 



$4 

$$ 

4 $ 




EI=1CC**SRI2/<VMAXUP+VS IK A) 
l*R I TE ( 6* 1 3G ) SRI2,E1 

130 FORMAT ( 1 8 X, 1 2HCRCSS RMS =,F9.4,9H FT/SEC =,F6.2»9H PER CENT! 

190 A1=E£AR/FC 

PRINT 15C,GAMMAC,GEFF,HC,HEFF f eBAR,Ai,VMAXUP 
150 F0RMAT(/1BX, »GAMMAC = ' f F6 .3 , 4X, * E FFEC 1 1 VE GAMMA =*,E6.3/18X% 

1 * HO/C = • * F6 . 3 ,4X » •EFFECTIVE H/D >* ,F6 . 3/1 8X f • BETA BAR *f, 

2F6.3.4X, 17HBETA * 0 = , F6 .3/38X, * MAX UPWASH VEL -MF6.2//1 

RETURN 
999 R MS-C • 

RETURN 

ENC 

SUER CUT I NE PRIM AT 

COMMON I YZC 0 C , 1 Y Z , NR # NC , A , B r C , Z C * V I N F f V 5 I N A , V Y ( 2 1 , 1 2 ) * V Z.( 21 » 1 2 ) * 
lY(3»12)*Z<211f RM S*NCRPF» I PR NT » NP T S » IRMS ,RM S V ZC ( 2 * 5 ) , IPRTM |NR AK 
CIMENSICN VC( 12 > , CH AfU 6 ) ,CHR ( 12 ) 

DATA CHAR/2H ,2H , 2F , 2H$$ , 1HY , 1HZ/ 

1 = 4+ IYZ 
CO 600 J= 1, NC 
600 VC ( J ) = 12 **Y ( 1 » J ) 

U = 1.6*RM5VZGU, 1 ) ♦ C * 5 

1FI IPRNT.EQ.2 > U = 0.1*VINF 

IFCI1.LT. 1) 11=1 

1 2=- II 

13=2*11 

A 1=11 

A2= 1 2 

A3= I 3 

WRITEI6, 610) CFARt I) t II, 12, I3,< VCl J) , J=l,NC > 

610 FORMAT!///, 49X, 37HMATRIX OF CROSS-SECT I CN VELGC l T IES//32XV48HCA 
1LCULATE0 VELCCITIES FRCP LEAST SQUARES FIT TO ,.A1, 22H COMPONENT 0 
2F VELOCITY///* Y,Z CENCTE Y,Z VELOCITY COMPONENT DATA*/ * C DENG 
3 TES THE CALCULATED VELCCITIES 1 //* INDICATES (V-VCALC) GREATER T 
4HAN*,I3,'. FT/SEC'/ • INDICATES (V-VCALC) LESS THAN* ,14, *. FT/S 




5EC* / * IS INDICATES AfiSI V-VC'ALC > GREATER THAN', 13,*. FT/ SEC* //IX ,13 
80UH*)/3X # 4HZ * , 51X , * AP PRGX . Y PCS IT ICNS */6X, 1H*,F 7. 2. I2F10.2) 

WR I T E C 6, 620 ) 

620 FCR R AT ( l X , 1 3C ( IE*) ) 

C DC LCCP 670 WRITES CUT THE VELOCITY FIELD DATA AND THE VELOCITIES USING FIT 
DO 67C 1 55 1 , NR 
ZW- 1 2 * *Z ( I ) 

K=( I+NCRPF + 6J/7 
IYZCCC=1 
CO 640 J = 1 , NC 
CHRC J)=CHAR( 1) 

VC< J) = F(Y4K, J),Z< W i 
IFCVYC I, JI.GT.1E3I GC TC 64 C 
El = V Y ( l » J )— VC ( J ) 

IF(ABS(E1)«LT.A1) GC TC 640 
I F ( A E S { E l ) • G T . A 3 ) GC TC 630 
IFIEl.GT.Al) C H R ( J ) = C H A R { 2 ) 

IFIEI.LT * A2 ) CHR < J ) = CHA R ( 3 ) 

GC TC 64C 

630 CHR ( J ) =CHAR ( 4 ) 

64C CONTINUE 

WRITE 16,660) <VYUtJ)»J = l*NCl 
WRITE16,6<30) ( VC ( J ) , CHR t J ) , J* 1 , NC ) 

I YZCCC=2 

CC 660 J- 1 , NC 

CHR ( J ) -CHAR (II 

VC( J)*FC YiK, J) ,ZU) ) —VS IN A 

IFCVZ.4I, J1.GT.1E3 J GC TC 66C 

El- VZ ( I, J }—VC( J ) 

I F ( ABSIEn.LT. All GC TC 660 
IF ( AES ( E 1 1.GT.A3 ) GC TC 65C 
IF(Ei.GT.Al) CHR ( J ) -CHAR ( 2 ) 

IFCE1.LT . A2 I CHR ( J I "CHAR ( 3 ) 

GC TC 660 






65C CHR(J)-CHAR(4) 

660 CONTINUE 

WRITE (6, ICO) ZMt(V2< l»J)f J=lfNC> 

670 WRITE(6*690) ( VC < J ) *CHR ( J ) , J- 1 * NC ) 

680 FCRHAT(6X,1H*/6X,1H*.12( IX,F6.1,*Y *)> 

690 FORMAT I6X* 1H** 12 1 1 X , F6» l * *C ' t A2 ) ) 

7C0 FGRMT(F6.2* 18+ • i 2 C IX,F6.1, *2 ' )) 

RETURN 

END 

FUNCTION F ( Y * Z Z ) , 

C FUNCTION F CALCULATES THE VELOCITY INCUCEC BY A PAIR OF DISTRIBUTED VORTICES. 

COHN ON IYZCOC>IYZ*NR»NC,A,B,C*ZO 

CON NON/DVC/F IB* F IC 

Z-ZZ-ZO 

E2-e*E 

YPC = Y + C 

YMC- Y~C 

R12=Z*Z+YPC*YPC 
R22=Z*Z+YMC*YRC 
IF (R12.E0.O. ) R12-1F— IC 
IF (R22.EC.0. I R22=1£-1C 
EXP1=EXP(-B2*R12) 

EXP 2- EXP (-B2*R22 ) 

El- ( I.-EXPl ) /R I 2 

E2=(1.-EXP2)/R22 

IF (RI2.LE.1E-1C) E 1 ^£2 

IF (R22.LE.-IE-1C) E2 = B2 

IF ( IYZCGC.EQ • 2 ) GO TO 1 

F-A*Z*(E2-Ell/6-283185 

RETURN 

1 F=A*(El*YPC~E2*YRC)/ 6.26316 5 
RETURN 

ENTRY OF ( Y * Z Z ) , 

C ENTRY CF CALCULATES THE PARTIAL DERIVATIVES LSEO IN A C I FFERENT I AL CORRECTION 



METhGC FOR A GAUSSIAN -DISTRIBUTED VORTEX PAIR. 

IF (R12.GT.1E— 1C ) DVhl=2^*YPC*< EXPi* ( i. + B2*R 12)-l. ) /HU2*K12) 

IF (R12.L6.1E-1C) DVP1=-2.*YPC*B2*B2 

IF (R22.GT.1E-1C) DVF2 = 2 .#YPC* ( EXP2*( 1 • + B2*R2 2 )-l . ) / IR22*R221 
IF (R22.LE.1E-1C) OVH2=-2.*YPC*B2*B2 
CF=F/ A 

IF ( IYZCCC. EG.2 ) GO TO 2 
FIB=A*Z*E*<EXP2-EXPl>/3. 141593 
F IC = -A*Z*(DVHl+CVF 2 )/fc. 2 8 3185 
RETURN 

FI8 = A*B* ( YPC*EXP l— YPC*E XP2 )/3. 141593 
F IC = A* (El+DVPl*YPC+E2+CVP2*YPC) /6. 2 8 3185 
RETURN 
ENC 

SUBROUTINE AUXIL 

PERFQRPS AUXILIARY OPERATIONS UPON INPUT DATA, 

COPPCN I YZCOD* IYZ,NR ,NC ♦ A » B t C » ZC* V I Nf » V S INA f VY ( 2 1 » 1 2 ) » VZ121 f 12 ) 
lYC3t 12)*Z(21 UHFS.NCRPFt IPRNT,NPTS f IRKS tRPSVZO ( 2 , 5 ) , IPRTK*NRAK 
CCPPCN/QNE/MKtNCATA, IRUN (4 ) ,R, AM JET , X4 t Z4 • SC , ANGVRT, VRANG 
C I PENS ION VC { 1 2 ) 

CATA AST/ 1H*/ 

IF ( IPRTM.EQ.C) GC TO 185 
WRITE (6, 150) 

15C FCRPAT (// 10X, * INPUT POSITION CATAG*) 

DC 170 1=1, NR 

K=(l4NCRPF+6)/7 

WRITE (6, 160) ( Y (K, J) , J= UNO > 

16C F 0 R P A T (■ / 1 2 ( F7.2, 1 Y 1 )) 

170 WR I TE ( 6 , 1 80 ) ( Z ( I) » J = 1 * NC > 

180 F OR P AT ( 1 2 ( F 7 . 2 » * Z • ) ) 

COUNT NUMBER OF USABLE DATA POINTS 
185 P=NR*NC 

CO 2C0 1=1, NR 

I F ( I • EG. 1 .OR. NRAK.EG.l) GC TO 190 




IF C A6S (Z < I)-Z( l- in .LT..5J M=M~NC 
190 CC 2C0 J*1,NC 

I F I V V I I»J)*LT-5CC.) GC TC 2CC 
F=M-1 

V Y(I,J)*IE6 
V Z 1 1 , J)*l-E6 
200 COM IMJE 

PRINT 2 1C t M 

2 1C FCRMT(/I4,» DATA POINTS') 

NC A T A-M 

C PERFORM VCRTICITY CALCULAT IONS 
A2-0 . 

I F ( IPRTM.EQ.l) PRINT 22 C 

220 FORMAT { / / * VCRTICITY CALCULATIONS'* 1CX* 24HUN I TSO UPPER - FT**2/SEC 
i /4G X , * LOWER - 1 /SEC * / ) 

CC 290 1*2, NR 

IF(IPRTM.EQ.l) PRINT 2 3C , < AST , J* 1 ,NC > 

230 FORM AT ( 5X, 12( A1 , 7X) ) 

Ai*2(I-lW( I) 

K*l I +NCRPF+6 ) /I 
L- ( I +NDRPF + 5 l/l 
CO 26C J *= 2, N C 

IF( VYU-ltJ-lt.GT.lE3) GC TC 240 
I F < V Y I I ,J-1)*GT.1E3> GC TC 240 
IF<VY(I-1.,J J.GT.1E3) GC TC 240 
I F ( V Y { I ,J ) . L T .1 E 3 ) GC TC 250 
240 VC ( J )* 1E6 
GC TC 260 

250 VC( J)*U VYf 1-1* J-U-VYf I , J ) >* ( Y < L , J )-Y ( K* J- 1 ) ) + { VY ( 1-1 , Jl-VYU , J-l 
1 ) >*<Y<K, J)-Y(L* J-l) ) + IVZ( [-1, J-1)+VZ( I, J-1)-VZ< 1-1* Jl-VZ H * J1 3 *AI) 
2/24 • 

A2=A2+VC(J) 

260 CCNTINUE 

IF < IPRTM *EQ* 1 ) PRINT 27C,(VC(J),J-2*NC) 


•*sj 

N3 





270 FCRMAT(5X*1241X*F7*3 ) ) 
rr pflG J = 2 * NC 

280 VCt J)=VC t J)*288./ABSI Al*(Y(K»J)+Y<Lt J)“V< K, J- l )- Y U , J-l > I * 

IF( IPRTM^EQ.ll PRINT 30C » ( VCI J ) * J = 2 »NC ) 

290 CONTINUE 

300 FORMAT (5X%12F8* 2 ) 

IRIPRTM.EQ.l) PRINT 23 C * ( AS T * J- l » NC ) 

A3 = 1 *5*A2/V INF 

PRINT 31C*A2*A3 _ #cc<~ r*UMA — cjl i i 

310 FORMAT ( / * SUMMED VORTICITY =»*F9.3*21H FT**2/SEC GAMMA -,F6.3I 

C DC LOOPS 320 AND 33C CONVERT POSITION CATA TO FEET 
DC 320 1=1* NR 
320 Z( l) = Z(I )/12. 

ZO=ZC / 12* 

CG 330 I=1,NRAK 
CC 3 30 J = 1 * NC 
330 Y(I*J)=Y(I*J)/12. 

RETURN 


C 


ENC 

SUBROUTINE QK A J { NOMR ,CM J * CMM * 0MN2 ) ypi rr * TV RATIO 

PURPOSEC CHCCSE THE MAJOR AND MINOR G CONDITIONS FOR THE GIVEN VELOCITY RATIO 

CM J ~C . 0 


GMM = 0 • 0 

GMN 2=0 *0 k . 

GO TC ( 28 * 28*21 * 22 * 23 * 24 * 25 * 26 * 28*27 )* NOMR 

21 CMJ=35.9 
GMN 1=15*5 
RETURN 

22 CM J = 2 C • 2 
CMM = 52.8 
CMN 2 = 8 »7 
RETURN 

23 C M J = 3 3 • 3 
RETURN 



24 CM J=23 • 1 
GHM = 35.0 
CMN2-9 *1 
RETURN 

25 CP J=26* l 
RETURN 

26 CM J - 1 9 • 9 
C PM -12 • 9 
CPN2=5.1 
RETURN 

27 C M J - 1 2 • 8 
CPM = 8 • 2 
CMN2 -3 • 2 

28 RETURN ^ 

ENC ^ 

FUNCTICN ARCLNC- ( X*R ) 

C PURPGSEC CALCULATE THE ARC LENGTH ALONG THE UNIFIED VORTEX PATH EQUATION* 

C GIVEN X/C AND R. 

CIMENSION G ( 8 } » V ( 8 ) 

DATA V/ .0950 125 ,.2816C35*. 4580167.. 6178762*. 7554 C44 
1 *.8656 31 !»• 9 44 5 7 50, .989 4009/ ,G/. 189 4506 ». 1826034 * 

2. 16 9 1565,. 1495959,. 1246289, . C95 1 585 « .0 6 2 2 53 5 « . 027 1524/ 

DATA E2» Al» A3/2. 660516* .463 3684 *-l . 141726/ 

E3= . 14i4484*R**2.253C?2*X**A3 
S*C. 

EG 1 I r 1 8 

1 S=S + G( I ) * (SGRTU l.+V ( I ) )**E2*E3) + SORT ( < 1 *-V ( I ) )**£2^E3J I 
ARC L NG= A 1*X* S 
RETURN 
ENC 

SUfiRCUT JIVE M I N V 3 ( A » C ) 

q PURPGSEC PERFORM PAT R IX INVERSE OF NORMAL EQUATION PAIR IX FROM LEAST SQUARES. 

C IRENS I C N A ( 6 ) » £ ( 3 ) 

REALMS A,6 f C*C 



E<l>*Al3)*A<5.)-mi*M4) 
B12)«AI3)*At Al-A(2)*A<5| 
e<t3)*A(2)*A(A)-A(3)*A(3) 

0~ A ( n*8Ui+A*2J*Bt 2 3) 

IFtC.EQ.C.) RETURN 

C =iA(l)*AI5)-M31M(3n/C 

A(5)=M AC 2 )*A«3)-A { )/l 

A (6) = < AU>*A<3)-A<2)*A< 2U/C 

A ( 1 ) = fi ( 1 F/D 

A ( 2 ) =6 ( 2 I /D 

A(3)-fi(3)/D 

A 1 4 ) “C 

RETURN 

ENC 
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